Este R Markdown recoge el enunciado de la práctica de la asignatura de redes sociales.

El objetivo es analizar un grafo, que se provee como fichero en el mismo paquete que este enunciado. En este fichero, encontramos solamente dos columnas, correspondiente a una interacción entre dos nodos de la red. Esta red está formada por distintos individuos que tienen contactos cara a cara durante un período de tiempo.

A continuación, dividimos la práctica en apartados, con una breve descripción de qué debe contener cada chunk de código donde el alumno desarrollará su respuesta así como las explicaciones que considere oportunas. Por favor, razona todas tus soluciones y escribe las explicaciones en azul.

Junto al título de cada apartado se encuentra la puntuación del mismo (pueden obtenerse hasta 10,5 puntos, aunque solamente se evaluará del 0 al 10).

Carga de datos y comprobaciones iniciales (0,5 puntos)

En este apartado, se pide:

library(igraph)
## Warning: package 'igraph' was built under R version 4.1.3
## 
## Attaching package: 'igraph'
## The following objects are masked from 'package:stats':
## 
##     decompose, spectrum
## The following object is masked from 'package:base':
## 
##     union
rm(list = ls())

setwd('C:/Users/Javier/Documents/masterAFI/18. Analisis de Grafos/01. Redes sociales/practica/')
# setwd('C:/Users/jherraez/Documents/masterAFI/18. Analisis de Grafos/01. Redes sociales/practica/')
dd <- read.csv('red_contactos.csv', sep = ';')

gg <- graph.data.frame(dd, directed = FALSE)
summary(gg)
## IGRAPH ad60b57 UN-- 1390 222744 -- 
## + attr: name (v/c)
vcount(gg)
## [1] 1390
ecount(gg)
## [1] 222744
gg2 <- simplify(gg, remove.multiple = TRUE, remove.loops = TRUE)

E(gg2)$weight = sapply(E(gg2), function(e) {
  length(all_shortest_paths(gg, from = ends(gg2, e)[1], to = ends(gg2, e)[2])$res) } )
summary(gg)
## IGRAPH ad60b57 UN-- 1390 222744 -- 
## + attr: name (v/c)
table(E(gg2)$weight)
## 
##     1     2     3     4     5     6     7     8     9    10    11    12    13 
## 12496 26438  2895  4818  1161  1481   642   739   364   416   278   246   173 
##    14    15    16    17    18    19    20    21    22    23    24    25    26 
##   193   134   140   106    97    94    71    73    65    54    65    45    45 
##    27    28    29    30    31    32    33    34    35    36    37    38    39 
##    40    30    37    27    38    23    23    35    18    11    18    19    20 
##    40    41    42    43    44    45    46    47    48    49    50    51    52 
##    11    15     9    13     9     7     9     7    11     5     8     9    10 
##    53    54    55    56    57    58    59    60    61    62    63    64    65 
##     7    12     5     8     7     7     3     7     7     5     6     6     7 
##    66    67    68    69    70    71    73    74    75    77    78    79    80 
##     5     4     4     2     3     3     1     1     1     3     1     1     1 
##    81    83    84    85    87    90    92    93    94    95    96    99   101 
##     1     1     2     3     2     1     3     1     1     1     1     2     1 
##   107   109   110   113   115   120   121   128   129   130   135   136   141 
##     1     1     1     1     1     1     1     1     1     1     1     1     1 
##   155 
##     1

Selección de la componente conexa mayor (0,5 puntos)

En este apartado, se pide realizar los pasos adecuados para generar un nuevo objeto grafo, que sea conexo, y que involucre a todos los nodos y enlaces de la componente conexa mayor del grafo original.

is.connected(gg2) # es conexo?
## [1] FALSE
ccs <- components(gg2)
ccs$csize
## [1] 1388    1    1
id_compmayor <- which.max(ccs$csize) # mayor componente conexa

vids <- ccs$membership == id_compmayor
vids <- which(ccs$membership == id_compmayor) #nodos de la comp conexa mayor

gg3 <- induced_subgraph(gg2, vids = vids ) #subbgrafo solamente con estos nodos

summary(gg3)
## IGRAPH e5c30ed UNW- 1388 53942 -- 
## + attr: name (v/c), weight (e/n)

Análisis descriptivo de la componente conexa mayor (2,5 puntos)

En este apartado, se pide analizar descriptivamente el grafo usando los conceptos que hemos visto durante las clases de teoría:

# grado medio
mean(degree(gg3))
## [1] 77.72622
hist(degree(gg3), breaks = 20)

# distancia media
average.path.length(gg3)
## [1] 3.066029
# diámetro
diameter(gg3)
## [1] 14
# Distribución de grados y ajuste a una Power-Law
graph.density(gg3)
## [1] 0.0560391
deg_dist <- degree_distribution(gg3)
fit <- fit_power_law(deg_dist)
fit
## $continuous
## [1] TRUE
## 
## $alpha
## [1] 3.447821
## 
## $xmin
## [1] 0.005043228
## 
## $logLik
## [1] 391.6632
## 
## $KS.stat
## [1] 0.1568648
## 
## $KS.p
## [1] 0.03535439
plot(deg_dist + 0.0001, log = 'xy', xlab = 'Node Degree', ylab = 'Probability')

lines(seq(deg_dist), seq(deg_dist)^-fit$alpha, col='red')

# Clustering
mean(count_triangles(gg3))
## [1] 1037.44
hist(count_triangles(gg3), breaks = 20)

# Entropía de los nodos Shannon
hist(diversity(gg3), breaks = 20)

# Centralidad de los nodos y comparación con métricas de grado y clustering
Degree <- degree(gg3)
Eig <- evcent(gg3)$vector
Closeness <- closeness(gg3)
Betweenness <- betweenness(gg3)

centralities <- cbind(Degree, Eig, Closeness, Betweenness)
round(cor(centralities), 2)
##             Degree  Eig Closeness Betweenness
## Degree        1.00 0.93      0.79        0.85
## Eig           0.93 1.00      0.69        0.78
## Closeness     0.79 0.69      1.00        0.59
## Betweenness   0.85 0.78      0.59        1.00

Análisis de comunidades de la componente conexa mayor (1,5 puntos)

En este apartado, se pide aplicar dos algoritmos de detección de comunidades, compararlos y seleccionar cuál es, en tu opinión, el que da una mejor respuesta. Razona tu selección.

#### infomap, fastgreedy, label propagation, walktrap

comms_infomap <- infomap.community(gg3)
comms_walktrap <- walktrap.community(gg3)

modularity(comms_infomap)
## [1] 0.4683851
modularity(comms_walktrap)
## [1] 0.4634219
compare(comms_infomap, comms_walktrap, method='nmi')
## [1] 0.812869
table(comms_infomap$membership)
## 
##   1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18  19  20 
##   9  55 392  45  26  11  37  25  12  30  61  87  22  52  57  60  35   4  13  18 
##  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40 
##  31  35   1  17  38  13   3  25   6  41  10  22  14   3  37  13   7   8   6   3 
##  41  42 
##   2   2
table(comms_walktrap$membership)
## 
##   1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18  19  20 
##  14  12  42  13 134 396  29  50  53  54  28  46   8  15  17  97   2   4  52  83 
##  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40 
##   5  10  16   7  30  13  19  32  17   2   2   1   1   1   1   1   1   1   1   1 
##  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  56  57  58  59  60 
##   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1 
##  61  62  63  64  65  66  67  68  69  70  71  72  73  74  75  76  77  78  79  80 
##   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1 
##  81  82  83  84  85  86  87  88  89  90  91  92  93  94  95  96  97  98  99 100 
##   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1 
## 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 
##   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1

Visualización del grafo por comunidades de la componente conexa mayor (1,5 puntos)

En este apartado, se pide visualizar el grafo coloreando cada nodo en función de la comunidad a la que pertenezca, según tu elección del apartado anterior.

wws <- ifelse(crossing(comms_infomap, gg3), 1, 500)
ll <- layout_with_fr(gg3, weights=wws)

palcol <- rainbow(n=max(comms_infomap$membership))

plot(gg3,
     layout = ll,
     vertex.label='',
     vertex.size=log(degree(gg3)),
     vertex.color=palcol[comms_infomap$membership])

Difundiendo un rumor (o un virus) en la componente conexa mayor (4 puntos)

Este apartado es el que más peso en la práctica tiene. Vamos a implementar un modelo epidemiológico sobre el grafo que, típicamente, se utiliza para simular escenarios de difusión de enfermedades pero también en contextos como la distribución de rumoeres e información. Vamos a implementar un modelo SIR que se caracteriza por tener los siguientes parámetros:

Se pide desarrollar una función que tenga como parámetros los tres valores anteriores y un cuarto que sea un grafo que, en nuestro caso, será la componente conexa mayor del grafo original de esta práctica. Dicha función simulará el proceso SIR:

Se pide ejecutar una simulación para tres o cuatro valores del parámetro beta (N y gamma pueden ser fijos en estas simulaciones) de este proceso de manera que se pueda visualizar:

generate_sir <- function(N = 5, beta = 0.1, gamma = 0.1, graph){
  graph <- set_vertex_attr(graph, 'label', value = 'S')
  
  sample_index <- sample(1:vcount(graph), N)
  graph <- set_vertex_attr(graph, 'label', index = sample_index, value = 'I')
  graph <- set_edge_attr(graph, 'label', value = '')
  
  n_contagiados = N
  n_susceptibles = vcount(graph) - N
  n_recuperados = 0
  dias_sin_contagios = 0
  dia = 1
  
  lista_contagiados = c(n_contagiados)
  lista_susceptibles = c(n_susceptibles)
  lista_recuperados = c(n_recuperados)
  lista_nuevos_contagios = c(n_contagiados)
  lista_grafo = list(graph)
  
  while (dias_sin_contagios < 3){
    n_contagiados_iniciales = n_contagiados
    n_contagiados_aux = n_contagiados
    n_nuevos_contagios = 0
    
    for (v in V(graph)){
      state <- V(graph)[v]$label
      
      if (state == 'I'){
        if (rbinom(n=1, size=1, prob = gamma) == 1){
          V(graph)[v]$label <- 'R'
          n_recuperados = n_recuperados + 1
          n_contagiados = n_contagiados - 1
        }
      }
      if (state == 'S'){
        for (neighbor in neighbors(graph, v)){
          if ((V(graph)[neighbor]$label == 'I') && (rbinom(n=1, size=1, prob = beta) == 1)){
            V(graph)[v]$label <- 'I'
            E(graph)[v %--% neighbor]$label <- 'I'
            n_contagiados = n_contagiados + 1
            n_contagiados_aux = n_contagiados_aux + 1
            n_susceptibles = n_susceptibles - 1
            break
          }
        }
      }
    }
    if (n_contagiados_iniciales == n_contagiados_aux){
      dias_sin_contagios = dias_sin_contagios + 1
    } else {
      dias_sin_contagios = 0
    }
    dia = dia + 1
    lista_contagiados = c(lista_contagiados, n_contagiados)
    lista_susceptibles = c(lista_susceptibles, n_susceptibles)
    lista_recuperados = c(lista_recuperados, n_recuperados)
    lista_nuevos_contagios = c(lista_nuevos_contagios, n_contagiados_aux - n_contagiados_iniciales)
    lista_grafo[[dia]] = graph
  }
  
  return (list('final_state' = V(graph)$label,
               'contagiados' = lista_contagiados,
               'susceptibles' = lista_susceptibles,
               'recuperados' = lista_recuperados,
               'nuevos_contagios' = lista_nuevos_contagios,
               'dias' = dia,
               'grafos' = lista_grafo
          ))

}

plot_evolution <- function(evolution){
  plot(1:evolution$dias, evolution$contagiados, type = 'l', col = 'red',
       ylab = 'Nº Personas', xlab = 'Días', ylim = c(0, max(evolution$susceptibles)))
  lines(1:evolution$dias, evolution$recuperados, type = 'l', col="green")
  lines(1:evolution$dias, evolution$susceptibles, type = 'l', col="blue")
  
  legend(evolution$dias, max(evolution$susceptibles)/2, 
         legend=c('Contagiados', 'Recuperados', 'Susceptibles'),
         col=c("red", 'green', "blue"), lty=1, cex=0.8, xjust = 1)
}

plot_infection_graph <- function(graph, index){
  
  graph_colors <- ifelse(V(graph)$label == 'I', 'red',
                       ifelse(V(graph)$label == 'R', 'green', 'blue'))

  subg <- subgraph.edges(graph, E(graph)[label == 'I'])
  
  plot(subg,
       layout = ll,
       vertex.label = '',
       vertex.size = log(degree(graph)),
       edge.label = '',
       vertex.color = graph_colors,
       main = paste('Dia', index))
  
}

create_gif <- function(list_graphs){
  img <- magick::image_graph(width = 1000, height = 1000)
  for (i in 1:length(list_graphs)){
    plot_infection_graph(list_graphs[[i]], i)
  }
  dev.off()
  animation <- magick::image_animate(img, fps = 2, optimize = TRUE)
  print(animation)
}

plot_new_infections <- function(list_new_infections, dias){
  plot(1:dias, list_new_infections, type = 'o', col = 'red', log = 'y', 
       ylab = 'Nº Personas', xlab = 'Días')
}
ejemplo_SIR_1 <- generate_sir(graph = gg3, beta = 0.1)

plot_evolution(ejemplo_SIR_1)

plot_new_infections(ejemplo_SIR_1$nuevos_contagios, ejemplo_SIR_1$dias)
## Warning in xy.coords(x, y, xlabel, ylabel, log): 5 y values <= 0 omitted from
## logarithmic plot

create_gif(ejemplo_SIR_1$grafos)
##    format width height colorspace matte filesize density
## 1     gif  1000   1000       sRGB  TRUE        0   72x72
## 2     gif   784    861       sRGB  TRUE        0   72x72
## 3     gif   613    855       sRGB  TRUE        0   72x72
## 4     gif   783    854       sRGB  TRUE        0   72x72
## 5     gif   621    854       sRGB  TRUE        0   72x72
## 6     gif   612    856       sRGB  TRUE        0   72x72
## 7     gif   614    855       sRGB  TRUE        0   72x72
## 8     gif   778    844       sRGB  TRUE        0   72x72
## 9     gif   778    845       sRGB  TRUE        0   72x72
## 10    gif   778    845       sRGB  TRUE        0   72x72
## 11    gif   777    845       sRGB  TRUE        0   72x72
## 12    gif   778    845       sRGB  TRUE        0   72x72
## 13    gif   738    852       sRGB  TRUE        0   72x72
## 14    gif   571    823       sRGB  TRUE        0   72x72
## 15    gif   779    855       sRGB  TRUE        0   72x72
## 16    gif   778    855       sRGB  TRUE        0   72x72
## 17    gif   438    819       sRGB  TRUE        0   72x72
## 18    gif   506    737       sRGB  TRUE        0   72x72
## 19    gif   734    777       sRGB  TRUE        0   72x72
## 20    gif   490    743       sRGB  TRUE        0   72x72
## 21    gif   524    718       sRGB  TRUE        0   72x72
## 22    gif   479    702       sRGB  TRUE        0   72x72

ejemplo_SIR_2 <- generate_sir(graph = gg3, beta = 0.05)

plot_evolution(ejemplo_SIR_2)

plot_new_infections(ejemplo_SIR_2$nuevos_contagios, ejemplo_SIR_2$dias)
## Warning in xy.coords(x, y, xlabel, ylabel, log): 10 y values <= 0 omitted from
## logarithmic plot

create_gif(ejemplo_SIR_2$grafos)
##    format width height colorspace matte filesize density
## 1     gif  1000   1000       sRGB  TRUE        0   72x72
## 2     gif   784    861       sRGB  TRUE        0   72x72
## 3     gif   596    861       sRGB  TRUE        0   72x72
## 4     gif   581    854       sRGB  TRUE        0   72x72
## 5     gif   591    854       sRGB  TRUE        0   72x72
## 6     gif   609    854       sRGB  TRUE        0   72x72
## 7     gif   611    855       sRGB  TRUE        0   72x72
## 8     gif   614    855       sRGB  TRUE        0   72x72
## 9     gif   612    846       sRGB  TRUE        0   72x72
## 10    gif   610    855       sRGB  TRUE        0   72x72
## 11    gif   610    844       sRGB  TRUE        0   72x72
## 12    gif   610    855       sRGB  TRUE        0   72x72
## 13    gif   610    788       sRGB  TRUE        0   72x72
## 14    gif   610    855       sRGB  TRUE        0   72x72
## 15    gif   551    757       sRGB  TRUE        0   72x72
## 16    gif   533    815       sRGB  TRUE        0   72x72
## 17    gif   544    739       sRGB  TRUE        0   72x72
## 18    gif   610    844       sRGB  TRUE        0   72x72
## 19    gif   578    789       sRGB  TRUE        0   72x72
## 20    gif   530    744       sRGB  TRUE        0   72x72
## 21    gif   610    844       sRGB  TRUE        0   72x72
## 22    gif   610    789       sRGB  TRUE        0   72x72
## 23    gif   610    861       sRGB  TRUE        0   72x72
## 24    gif   534    715       sRGB  TRUE        0   72x72
## 25    gif   565    783       sRGB  TRUE        0   72x72
## 26    gif   322    700       sRGB  TRUE        0   72x72
## 27    gif   542    722       sRGB  TRUE        0   72x72
## 28    gif   270    722       sRGB  TRUE        0   72x72
## 29    gif   581    789       sRGB  TRUE        0   72x72
## 30    gif   382    815       sRGB  TRUE        0   72x72
## 31    gif   293    714       sRGB  TRUE        0   72x72
## 32    gif   562    783       sRGB  TRUE        0   72x72
## 33    gif   301    715       sRGB  TRUE        0   72x72
## 34    gif   324    583       sRGB  TRUE        0   72x72
## 35    gif   294    705       sRGB  TRUE        0   72x72

ejemplo_SIR_3 <- generate_sir(graph = gg3, beta = 0.01)

plot_evolution(ejemplo_SIR_3)

plot_new_infections(ejemplo_SIR_3$nuevos_contagios, ejemplo_SIR_3$dias)
## Warning in xy.coords(x, y, xlabel, ylabel, log): 7 y values <= 0 omitted from
## logarithmic plot

create_gif(ejemplo_SIR_3$grafos)
##    format width height colorspace matte filesize density
## 1     gif  1000   1000       sRGB  TRUE        0   72x72
## 2     gif   784    861       sRGB  TRUE        0   72x72
## 3     gif   539    743       sRGB  TRUE        0   72x72
## 4     gif   556    798       sRGB  TRUE        0   72x72
## 5     gif   562    835       sRGB  TRUE        0   72x72
## 6     gif   567    806       sRGB  TRUE        0   72x72
## 7     gif   575    823       sRGB  TRUE        0   72x72
## 8     gif   571    855       sRGB  TRUE        0   72x72
## 9     gif   580    855       sRGB  TRUE        0   72x72
## 10    gif   578    855       sRGB  TRUE        0   72x72
## 11    gif   575    855       sRGB  TRUE        0   72x72
## 12    gif   588    854       sRGB  TRUE        0   72x72
## 13    gif   569    855       sRGB  TRUE        0   72x72
## 14    gif   579    854       sRGB  TRUE        0   72x72
## 15    gif   577    854       sRGB  TRUE        0   72x72
## 16    gif   574    854       sRGB  TRUE        0   72x72
## 17    gif   578    854       sRGB  TRUE        0   72x72
## 18    gif   579    854       sRGB  TRUE        0   72x72
## 19    gif   573    855       sRGB  TRUE        0   72x72
## 20    gif   588    855       sRGB  TRUE        0   72x72
## 21    gif   575    855       sRGB  TRUE        0   72x72
## 22    gif   549    746       sRGB  TRUE        0   72x72
## 23    gif   583    854       sRGB  TRUE        0   72x72
## 24    gif   546    765       sRGB  TRUE        0   72x72
## 25    gif   579    856       sRGB  TRUE        0   72x72
## 26    gif   573    856       sRGB  TRUE        0   72x72
## 27    gif   373    735       sRGB  TRUE        0   72x72
## 28    gif   562    799       sRGB  TRUE        0   72x72
## 29    gif   525    706       sRGB  TRUE        0   72x72
## 30    gif   420    701       sRGB  TRUE        0   72x72
## 31    gif   569    855       sRGB  TRUE        0   72x72
## 32    gif   175    605       sRGB  TRUE        0   72x72
## 33    gif   568    855       sRGB  TRUE        0   72x72
## 34    gif   195    579       sRGB  TRUE        0   72x72
## 35    gif   384    488       sRGB  TRUE        0   72x72
## 36    gif   316    738       sRGB  TRUE        0   72x72